Probing Tissue Multifractality Using Wavelet based 

Multifractal Detrended Fluctuation Analysis: 

Applications in Precancer Detection 
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Abstract — The refractive index fluctuations in the connective 
tissue layer (stroma) of human cervical tissues having different 
grades of precancers (dysplasia) was quantified using a wavelet- 
based multifractal detrended fluctuation analysis model. The 
results show clear signature of multi-scale self-similarity in the 
index fluctuations of the tissues. Importantly, the refractive index 
fluctuations were found to be more anti-correlated at higher 
grades of precancers. Moreover, the strength of multifractality 
was also observed to be considerably weaker in higher grades of 
precancers. These results were further complemented by Fourier 
domain analysis of the spectral fluctuations. 

Index Terms — Fractal, Multifractal detrended fluctuation anal- 
ysis, Phase contrast imaging, Refractive index, Cancer diagnosis. 



I. Introduction 

Over the last fifty years, there have been tremendous ad- 
vances in our understanding of the molecular and cellular 
processes of cancer; however, it still remains to be the deadliest 
disease of our time. Despite the significant progress made in 
treatment of a number of neoplastic disorders, early detection 
of neoplastic changes appears to be our best method to 
improve patient quality of life and reduce cancer mortality. 
The conventional methods for early detection and diagnosis 
of cancer rely on histological and cytological examination 
of tissue. This approach rely on characterizing morphologic 
and architectural alterations, including increased nuclear size, 
increased nucleus/cytoplasm ratio, hyperchromasia, pleomor- 
phism, and loss of normal epithelial architecture. There is 
a tremendous need in developing better methodologies for 
extraction and quantification of the morphological alterations 
associated with cancer/pre-cancer development. In this regard, 
the optical spectroscopic and imaging approaches have shown 
early promise in quantifying both the morphological (using 
elastic scattering spectroscopy) and biochemical (using in- 
elastic scattering spectroscopy such as fluorescence and Ra- 
man) alterations associated with cancer development [l]-[3]. 



Notably, polarized elastic scattering spectroscopy has been 
explored to quantify the self-similar (fractal) nature of micro- 
scale fluctuation of local refractive index in tissues. These 
studies have revealed that changes in tissue self-affinity can 
serve as a potential bio-marker for pre-cancer [4]. 

Since most of the cancers arise in epithelial tissues, majority 
of the previous attempts on developing methods (optical or 
non-optical) for early diagnosis of cancer relied largely on 
quantifying the alterations in the superficial epithelial layer. 
However, recently it has been recognized that in addition to the 
alterations in the superficial epithelial cells, neoplasia is also 
associated with characteristic changes in the underlying con- 
nective tissue layer (stroma). Progression of cancer involves 
complex interactions between neoplastic cells and the stroma. 
Also, carcinogenesis results, in part, from defective epithelial- 
stromal communication [5]. In fact, alterations in stromal 
biology may precede and stimulate neoplastic progression in 
pre-invasive disease [6]. Interestingly, the collagen fiber net- 
work present in stroma also exhibits fractal architecture in the 
organization of the fibers and micro-fibrils [7]. Quantification 
of the changes in the fractal characteristics of the stroma may 
thus provide additional targets to aid in screening and early de- 
tection of precancerous changes. With this motivation, we have 
explored the use of a wavelet-based multifractal detrended 
fluctuation analysis model to extract and quantify the fractal 
properties of the refractive-index structure (inhomogeneities) 
in the stromal layer of dysplastic (pre-cancerous) human 
cervical tissues. The results showed interesting variations in 
the micro-optical tissue fractal properties in the various states 
of dysplasia. The details of these results are presented in this 
paper. 

II. EXPERIMENTAL METHODS AND MATERIALS 

We have used a phase contrast microscope to measure the 
spatial fluctuations of the refractive index structures of the 



tissues. The samples used in this study were pathologically 
characterized, unstained sections of human cervical tissues 
(containing both epithelial and stromal regions) in glass 
slide. Tissue specimens were frozen and vertically sectioned 
(epithelium to stromal region) to a thickness of 5/j,m for 
measurements and analysis. The lateral dimensions of tissue 
sample were typically 4mm x 6mm. The dysplastic tissues 
were histopathologically characterized as Grade I, II and III 
respectively, based on their epithelial architecture. 

The phase contrast images were recorded from the stromal 
region of the tissue sections using a phase contrast microscope 
(Olympus 1X81, USA). Images of the specimens were taken 
at a magnification of 100X and were recorded with a CCD 
camera (ORCA-ERG, Hamamatsu) having 1344 x 1024 pixels 
(pixel dimension 6.45/im). 

III. Theory 

Biological tissue is an optically inhomogeneous medium 
having random fluctuations in local refractive indices (arising 
from the presence of microscopic inhomogeneities; macro- 
molecules, cell organelles, organized cell structure, extra- 
cellular matrix, interstitial layers etc.), the spatial scale of 
which varies from few dozen nanometers to several tens 
of micrometers. Tissue can thus be described as a random 
continuum of the inhomogeneities of the refractive index with 
varying spatial scale. As noted previously, for many types of 
tissues, the spatial distribution of refractive index have been 
found to exhibit statistical self-similarity [4], [8]. Specifically, 
the stromal tissues (investigated in this study) are known to 
be comprised of complex fibrous network having fractal-like 
arrays of collagen molecule, micro-fibrils and fiber bundles 
[7]. Since, the phase contrast image is a direct indicator of such 
spatial fluctuation of local refractive index structure, analysis 
of this fluctuation through appropriate statistical models may 
yield wealth of interesting micro-optical parameters, each of 
which can potentially serve as a useful biological metric. 
In order to explore this possibility, we have first unfolded 
(in one direction) the recorded images to obtain statistically 
large enough refractive index fluctuation series for analysis. 
These were then analyzed through (i) Fourier analysis and (ii) 
Wavelet based Multi Fractal De-trended Fluctuation Analysis 
(WB-MFDFA). 

A. Fourier Domain Analysis 

Fourier domain analysis is a simple and convenient ap- 
proach to test and quantify the fractal behavior of any fluc- 
tuation series. Such analysis based on the power spectrum of 
the Fourier transform of fluctuations has been widely used in 
many areas dealing with noise, from chaos to order transitions. 
Recently this approach has also been used to analyze light 
scattering from biological tissues [4], [8]. The power spectrum 
can be related to the Hurst scaling exponent (H, a measure 
of the fractal nature: < H < 1 for fractals, ) of any 
(one dimensional) fluctuation series exhibiting statistical self- 



similarity by the following relationship [9] 

P(k) « k~ a 

a = 2H + l 



(1) 



Here, in the power spectrum P(k), k is the frequency (spatial 
frequency in our case). 

Note that Eq. 1 is based on the assumption of an exact 
power law autocorrelation function of the fluctuations. In 
fact, such dependence can be understood by noting that the 
power spectrum of any fluctuation series is related to the 
autocorrelation function through an inverse Fourier transform. 
The exact form of the power spectrum is different though 
for fluctuations having more general type of autocorrelation 
functions (such as the von Karman self-affine function, which 
has been extensively used to describe refractive index corre- 
lations in tissues and other similar type of random media [4], 
[8]). Never-the-less at the limit of large spatial frequencies, 
the power spectrum for those also essentially converge to Eq. 
1[8]. 

In the first step, we have thus analyzed the recorded tissue 
refractive index fluctuations (one dimensional) employing dis- 
crete Fourier transform (DFT) and its power spectrum. The 
obtained power spectrum was fitted to Eq. 1 to yield the value 
for the exponent a and the Hurst scaling exponent H. 

Note that the above analysis is based on monofractal hy- 
pothesis, which assumes that the scaling properties are the 
same for the entire region of the tissue section examined. 
However, considering the wide range of the dimensions of 
the inhomogeneities (in refractive index) and the complex 
nature of the correlations present between them in tissue (the 
stromal tissues), this assumption may be unrealistic. In such 
situation, it would be desirable to use a more general type 
of analysis which could extract and quantify the nature of 
multi-fractality (if there) present in the signal. A multifractal 
signal can be decomposed into many subsets, wherein a 
local Hurst exponent quantifies the local singular behavior 
and gives the local scaling of the signals. Such multifractal 
analysis may provide additional diagnostic information on the 
nature of the nonlinearity encoded in the Fourier space of 
the refractive index fluctuations. We therefore applied a more 
general wavelet based approach on the multifractal detrended 
fluctuation analysis (WB-MFDFA), as discussed below. 

B. Wavlet Based Multi Fractal Detrended Fluctuation Analy- 
sis (WB-MFDFA) 

The multi fractal detrended fluctuation analysis (MFD FA) 
has been introduced to study the long range correlations in 
fluctuations when a trend is present [10]. It has been found to 
be very useful for studying non-stationary fluctuation series, 
where it is of vital importance to separate local fluctuations 
from average behavior (trend). In MFDFA and its variants, the 
detrending is usually done using polynomial fits. However, a 
more recent approach exploits the use of discrete wavelets 
for performing the detrending procedure [11], [12]. It has 
been shown that the natural, built-in variable window size 
in wavelet transforms makes this procedure more efficient 



for this purpose. We have therefore explored the use of this 
Wavelet based Multi Fractal Detrended Fluctuation Analysis 
(WB-MFDFA) for analyzing our phase contrast data from 
tissue sections [11]. 

In this approach, the profile is first calculated from a series 
Xt(t = 1, . . . , N) of length N, by performing the cumulative 
sum of the series after subtracting the mean. 



Y(i) = 5> t - <s>], i = l,...,N. 



(2) 



In order to separate the fluctuation from the trend, a discrete 
wavelet transform is then carried out on the profile Y(i). 
Discrete wavelets belonging to Daubechies family (Db4) is 
used for this purpose. It may be worth mentioning here that 
since, these wavelets satisfy the vanishing moment conditions, 
the low-pass coefficients keep track of the polynomial trend in 
the data. Thus, reconstruction using the low-pass coefficients 
alone is quite accurate in extracting the local trend, in a 
desired window size. The fluctuations are then extracted at 
each level by subtracting the obtained series from the original 
data. The details of this detrending procedure can be found in 
[11]. Thus obtained fluctuations are then considered for further 
analysis, as follows. The extracted fluctuations are subdivided 
into non-overlapping segments M s = int(iV/s) where the 
window size at a particular level (L) is, s = 2 i ~ L ~ 1 'W for 
the chosen wavelet. W is the number of filter coefficients 
of the discrete wavelet transform under consideration. This 
procedure is repeated from the end to beginning to calculate 
the local variance in case N/s is not an integer which would 
mean that some data points have to be discarded. Following 
this procedure, the q t h order fluctuation function, F q (s) is 
obtained as, 

( 2M 3 \ V8 

-E^ 2 (Mf /2 
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Here 'q' is the order of the moment which can take any real 
value. This procedure is repeated for various values of q (q ^ 
0) for variable window sizes (in our case we varied — 10 < q < 
10). The scaling behavior of the fluctuations can be obtained 
by analyzing the fluctuation function in a logarithmic scale for 
each value of q 



F q (s) 



Mq) 



(4) 



Note that special care must be taken for a value of q = 0. 
This follows because for q = 0, direct evaluation of Eq. 
3 leads to divergence of the scaling exponent [11]. In that 
case, logarithmic averaging has to be employed to find the 
fluctuation function 
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(5) 



The generalized Hurst exponent h(q) values are independent 
of q if the series is mono-fractal. For multifractal series, h(q) 
values depend on q. Note that the value for h at q = 2, 
corresponds to the Hurst exponent (H = h(q = 2)) for a 



stationary series. In fact, the correlation behavior of any series 
(stationary or non-stationary) is characterized from the Hurst 
exponent. For long range correlation, H > 0.5, H = 0.5 
for uncorrected and H < 0.5 for long range anti-correlated 
[9]. The variation h(q) can thus be used as an indicator to 
correlation behavior and multifractality of the fluctuations. 

IV. Results and discussions 

Typical phase contrast images recorded from the stromal 
region of dysplastic cervical tissues, histopathologically char- 
acterized Grade I and Grade III, are shown in Figure la and 
lb respectively. The wide range of sizes and shapes of the 
index inhomogeneities and their high packing densities are 
evident. The complex nature of the index variation however, 
underscores the problem how does one extract and quantify 
useful diagnostic metrics from this intertwined information? 
The observed index variations were then unfolded to generate 
the fluctuation series, which were subsequently analyzed via 
the Fourier analysis and WB-MFDFA. 

A. Results of the Fourier analysis 

Figure 2(a) and 2(b) show the corresponding power spectral 
density of the index inhomogeneities in the two tissues, which 
we obtained by performing Fourier transformation on the 
unfolded fluctuation series. The power spectral density and 
the spatial frequency are displayed in a log-log scale in these 
figures. Several interesting characteristics of the spectrum can 
be discerned. Beyond a certain spatial frequency range, the 
spectral density appears linear when plotted on a loglog scale, 
indicating a power-law scaling which depends on the sizes of 
the inhomogeneities that contribute to index variations. How- 
ever, this is also associated with large background fluctuations 
indicating the overall randomness of the underlying index 
variations. Interestingly, the power law coefficient does not 
appear to be uniform throughout the entire spatial frequency 
region. This behavior is a clear indication of the multifractal 
nature of the index fluctuations. The data were therefore fitted 
to Eq. 1 for two different selected spatial frequency range; 
(i) for the overall broad range of the spatial frequencies where 
the linear behavior was noted (power law exponent noted as a 
in the figures) and (ii) at a narrower higher spatial frequency 
range where the slope was apparently different (power law 
exponent noted as a\ in the figures). The results of these 
fits are displayed as solid and dashed lines respectively, in 
the figures. The corresponding estimates for the exponents 
are also listed in the figures. Importantly, the value for the 
estimates for the average Hurst parameter are significantly 
different, H = 0.416 and 0.224 for Grade I and Grade III 
dysplastic connective tissues respectively. A lower value of 
H is indicative of increasing roughness of the medium, i,e 
predominance of index inhomogeneities having smaller spatial 
dimensions. Note that the limiting values of H, unity and 
zero correspond to a smooth Euclidean random field (marginal 
fractal) and a space-filling field (extreme fractal) respectively 
[9]. The lower value of H in higher grades of dysplasia is 
consistent with previous reports where it has been observed 
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Fig. 1. (Color Online) Phase-contrast images of the stromal region of cervical tissues 
having (a) Grade I and (b) Grade III stage of dysplasia. 
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Fig. 2. (Color Online) Power spectra of index variations in the stromal region of 
cervical tissues having (a) Grade I and (b) Grade III stage of dysplasia. The power 
spectral density and the spatial frequency are displayed in a log-log scale. 



that the volume fraction of collagen fibers in the stroma 
decreases with precancer and cancer progression, and fibers 
tend to be shorter and more disconnected in neoplastic stroma 
[7]. 

Never-the-less the observed nonlinearity in the Fourier space 
of the refractive index fluctuations (different values of the 
exponent a at different spatial frequency range) underscores 
the need for employing a multi-resolution fluctuation analysis 
model. 

B. Results of the analysis of WB-MFDFA 

In Figure 3(a) and 3(b), we summarize the results of the 
Wavelet based Multi Fractal Detrended Fluctuation Analysis 
performed on the unfolded index fluctuations series in the two 
tissues (whose results have been presented in Figure 1 and 2). 
The dependence of the generalized Hurst exponent h(q) on 
the order of the moment (q) are displayed in the figures. 

As is apparent, the index fluctuations for both these tis- 



sues exhibit long range correlations and multifractal behavior 
(manifested as a variation of h(q)). Interestingly, the variation 
of h(q) is observed to be considerably weaker for the Grade 
III dysplastic connective tissue, as compared to that for the 
Grade I tissue. This clearly indicates that with increasing 
grades of dysplasia (with precancer progression), the strength 
of multifractality also decreases. This was further quantified 
and confirmed by determining the singularity spectrum of 
multifractality (not shown here) [10], [11]. As expected, the 
width of the singularity spectrum was found to be considerably 
narrower in case of the Grade III tissues compared to the Grade 
I tissues, indicating a trend towards monofractal behavior. 
The values for h(q = 2) for these two tissues (Grade I and 
III) were determined to be 0.55 and 0.203 respectively. Note 
that these values should correspond to the Hurst exponent 
(H = h(q = 2)) for a stationary series. The observed 
variations in the values for these parameters (in the two 



tissues) are qualitatively similar to that observed for the Hurst 
exponent (H) determined using the Fourier analysis (in both 
cases, the value for the Hurst exponent decreases at higher 
grade of precancer indicating a trend towards long range anti- 
correlation). The absolute values estimated using these two 
approaches are different. This is expected because, the Fourier 
analysis estimates this parameter (H) using a monofractal 
approximation that is far from being true for these tissues (as 
is apparent from the non-stationary nature of the fluctuations). 
Importantly, none of the multifractal trends observed in these 
tissues could be gleaned using the Fourier analysis, where 
at best one can qualitatively predict about the overall nature 
of the self-similar behavior of the tissue refractive index 
fluctuations. Derivation and quantification of the multifractal 
trends is possible due to the multi-resolution ability of the 
WB-MFDFA approach. 
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Fig. 3. (Color Online) The variation of generalized Hurst exponent h{q) as a function 
of order of the moment (q) for the index variations in the stromal region of cervical 
tissues having (a) Grade I and (b) Grade III stage of dysplasia. The value for h(q — 2) 
for these two tissues were 0.55 and 0.203 respectively. 



V. CONCLUSIONS 

To conclude, we have employed a wavelet based ap- 
proach on the multifractal detrended fluctuation analysis (WB- 
MFDFA) to extract and quantify self-similarity of tissue 
refractive index fluctuations. This novel approach for fluc- 
tuation analysis has been explored here in a scenario of 
significant clinical interest, that for detecting precancerous 
alterations in human cervical tissues. Early indications show 
promise of this approach for extraction and quantification of 
the morphological alterations associated with precancers. The 
information obtained on the multifractal nature of refractive 
index fluctuations in tissue may also prove to be valuable for 
developing light scattering based spectroscopic and imaging 
techniques for non-invasive and early diagnosis of cancer. This 
is currently under investigation in our laboratory. 
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